Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Support write and read of boolean masks in wcs1d-fits #1051

Merged
merged 3 commits into from
Apr 18, 2023

Conversation

dhomeier
Copy link
Contributor

This PR includes a test that was missed before merging #1009 and adds functionality to write Spectrum1D masks of dtype = bool in WCS1D format, resolving an issue that came up in #1009 (comment).
As bool or bit cannot be written natively to FITS, enabling round-tripping requires encoding the original type, for which there also seems to be no standard. I propose to introduce the header keyword 'BFORM', which is neither defined nor reserved in the FITS standard, and use it analogously to 'TFORMn' in BINTABLE extensions, for which these values are allowed (excerpt):

TFORMn value Description
’L’ Logical
’X’ Bit
’B’ Unsigned byte
’I’ 16-bit integer

Here BFORM is only used for the 1-bit types, as other formats are already recognised, but it might be added for others for consistency.

@dhomeier dhomeier added the io label Mar 29, 2023
@dhomeier dhomeier requested a review from pllim March 29, 2023 15:22
@dhomeier dhomeier added Extra CI Run cron in PR and removed Extra CI Run cron in PR labels Mar 29, 2023
@dhomeier
Copy link
Contributor Author

@rosteen do PRs need extra approval for the CI?

@rosteen
Copy link
Contributor

rosteen commented Mar 29, 2023

@rosteen do PRs need extra approval for the CI?

They shouldn't, since you've contributed here before.

Copy link
Member

@pllim pllim left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM. Thanks!

@pllim
Copy link
Member

pllim commented Mar 29, 2023

Actions was having problem. I am going to try to close/re-open.

@pllim pllim closed this Mar 29, 2023
@pllim pllim reopened this Mar 29, 2023
@pllim
Copy link
Member

pllim commented Mar 29, 2023

If CI does not kick off still, maybe just wait out the server problems. https://www.githubstatus.com/

@dhomeier
Copy link
Contributor Author

They shouldn't, since you've contributed here before.

And I could even trigger Extra CI :) but workflows still seem to be clogged...

@codecov
Copy link

codecov bot commented Mar 29, 2023

Codecov Report

Attention: Patch coverage is 88.88889% with 2 lines in your changes missing coverage. Please review.

Please upload report for BASE (main@2a7e5aa). Learn more about missing BASE report.
Report is 80 commits behind head on main.

Files with missing lines Patch % Lines
specutils/io/default_loaders/wcs_fits.py 88.88% 2 Missing ⚠️
Additional details and impacted files
@@           Coverage Diff           @@
##             main    #1051   +/-   ##
=======================================
  Coverage        ?   70.29%           
=======================================
  Files           ?       64           
  Lines           ?     4410           
  Branches        ?        0           
=======================================
  Hits            ?     3100           
  Misses          ?     1310           
  Partials        ?        0           

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

pllim added a commit to pllim/jdaviz that referenced this pull request Mar 30, 2023
to write out fitted cube.

Need astropy/specutils#1051

[ci skip] [rtd skip]
@pllim
Copy link
Member

pllim commented Mar 30, 2023

Not quite there. Looks like flux is written to EXT 0 even when there is mask, which is super confusing. Usually when there are multiple data extensions, EXT 0 is primary. Let's see this example:

import numpy as np
from astropy import units as u
from astropy.io import fits
from astropy.wcs import WCS

from specutils import Spectrum1D

w = WCS(naxis=3)
w.wcs.ctype = 'FREQ', 'DEC--TAN', 'RA---TAN'
w.wcs.set()

sp = Spectrum1D(flux=np.ones((8, 9, 10)) * u.nJy, mask=np.zeros((8, 9, 10), dtype=np.int8), wcs=w)
sp.write("ztmp.fits", format="wcs1d-fits", overwrite=True)
fits.info("ztmp.fits")

With this patch, I get this:

No.    Name      Ver    Type      Cards   Dimensions   Format
  0  FLUX          1 PrimaryHDU      29   (10, 9, 8)   float64
  1  MASK          1 ImageHDU        11   (10, 9, 8)   uint8 (rescales to int8)

But I want this:

No.    Name      Ver    Type      Cards   Dimensions   Format
  0  PRIMARY       1 PrimaryHDU     303   ()
  1  FLUX          1 ImageHDU        29   (10, 9, 8)   float64
  2  MASK          1 ImageHDU        11   (10, 9, 8)   uint8 (rescales to int8)

pllim added a commit to pllim/jdaviz that referenced this pull request Mar 30, 2023
to write out fitted cube.

Need astropy/specutils#1051

[ci skip] [rtd skip]
@dhomeier
Copy link
Contributor Author

Not quite there. Looks like flux is written to EXT 0 even when there is mask, which is super confusing. Usually when there are multiple data extensions, EXT 0 is primary. Let's see this example:

Yes, if you want the flux written to a different extension, you need to specify (e.g.) hdu=1.
The default has always been to write flux to the first available HDU, which for ImageHDU format is the primary HDU.
Mask and uncertainty, if present, are written to the subsequent HDUs per the recent PRs, but flux will be in the same HDU it has always been. I don't quite see how this becomes more confusing when a mask HDU is appended – flux is in HDU 0, whether mask and/or uncertainty are there or not.

@pllim
Copy link
Member

pllim commented Mar 31, 2023

flux is in HDU 0, whether mask and/or uncertainty are there or not

But this is not how 90% of the data out there are laid out. I have never seen a real life example where the FLUX is still in primary extension if there are ERR or DQ tied to it.

@dhomeier
Copy link
Contributor Author

But this is not how 90% of the data out there are laid out. I have never seen a real life example where the FLUX is still in primary extension if there are ERR or DQ tied to it.

The SDSS spPlate format e.g., with FLUX in the primary HDU and InvVar in the first extension.
SDSS-I/II spSpec also has FLUX in the primary HDU, though with ERR and MASK etc. folded into the same data "cube", and the generic IRAF reader does not even support reading from anything but HDU 0 (but also does not account for any uncertainty or mask data).
APOGEE uses a scheme more similar to what you are describing, but writing flux, StdErr and mask to extensions 1, 2 and 3, but it has its own dedicated loader – however perhaps wcs1d should adhere to that order as well, writing uncertainty before mask?
Bottom line, I don't have any robust percentage numbers on how all spectra out there are distributed, but for the basic image format using the primary HDU seems at least as common as the first extension. In any case the default for writing a flux-only spectrum has been hdu=0 since specutils 1.1, and I don't think it's a good idea to change that now.
OTOH changing the default HDU to 1 depending on the presence of uncertainty and/or mask seems rather more confusing to me (e.g. what should it be if a mask exists, but is skipped in writing by mask_name=None?). But if that is what the majority of users find more logical, so be it...

@pllim
Copy link
Member

pllim commented Mar 31, 2023

What do the other specutils maintainers think?

should adhere to that order as well, writing uncertainty before mask?

Yes, I think it is common to have the SCI, ERR, DQ in that order.

@dhomeier
Copy link
Contributor Author

dhomeier commented Apr 4, 2023

Yes, I think it is common to have the SCI, ERR, DQ in that order.

I've changed the order; when updating the tests I realised that the current wcs1d-fits auto-detection also expects the science data in the primary HDU, as it is only probing data structures in HDU 0. So changing this would also break something – either identifying previously written files, or new files (with mask)... In principle identify_wcs1d_fits could also be extended to search through all HDUs for a match, but that would introduce additional complexity (and possibly also performance penalties).

@rosteen rosteen requested review from keflavich, eteq, nmearl and rosteen April 5, 2023 21:02
pllim added a commit to pllim/jdaviz that referenced this pull request Apr 7, 2023
to write out fitted cube.

Need astropy/specutils#1051

[ci skip] [rtd skip]
Copy link
Member

@pllim pllim left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I had to pass in extra options (hdu=1, flux_name="SCI") to make an existing test pass over at spacetelescope/jdaviz#2094 but it is not the end of the world.

LGTM. Thanks!

Copy link
Contributor

@rosteen rosteen left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Code looks fine to me, thanks!

@rosteen rosteen merged commit b332b87 into astropy:main Apr 18, 2023
@pllim
Copy link
Member

pllim commented Apr 18, 2023

Please give me a ping when this ends up in a released version, so I can follow up over at Jdaviz to use this for Cubeviz. Thanks!

pllim added a commit to pllim/jdaviz that referenced this pull request May 10, 2023
to write out fitted cube.

This uses astropy/specutils#1051

Update specutils pin.
@pllim
Copy link
Member

pllim commented May 10, 2023

When is the next release for specutils?

@dhomeier dhomeier deleted the wcs1d-masked branch May 10, 2023 18:34
rosteen pushed a commit to rosteen/specutils that referenced this pull request Aug 15, 2023
* Test warning on unspecified uncertainty type

* Store and retrieve boolean masks in wcs1d-fits

* Move uncertainty before mask in hdulist
pllim added a commit to pllim/jdaviz that referenced this pull request Jul 29, 2024
pllim added a commit to pllim/jdaviz that referenced this pull request Jul 30, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants